#include <iostream>
#include <vector>
#include <cmath>
#include <fstream>

// 定义原始信号函数
double originalSignal(double t) 
{
    // 正弦波作为示例信号
    return std::sin(2 * M_PI * 10 * t);
}

// 采样函数
std::vector<double> sampleSignal(double (*signal)(double), double samplingRate, double duration) 
{
    std::vector<double> samples;
    double dt = 1.0 / samplingRate;
    for (double t = 0; t < duration; t += dt) 
    {
        samples.push_back(signal(t));
    }
    return samples;
}

// 重建信号
double reconstructSignal(const std::vector<double>& samples, double samplingRate, double t) 
{
    double dt = 1.0 / samplingRate;
    int index = static_cast<int>(t / dt);
    if (index < samples.size()) 
    {
        return samples[index];
    }
    return 0.0;
}

// 计算重建信号与原始信号的均方误差
double calculateMSE(double (*original)(double), const std::vector<double>& samples, double samplingRate, double duration) 
{
    double mse = 0.0;
    double dt = 1.0 / 1000;  // 采样间隔
    int numPoints = static_cast<int>(duration / dt);
    for (int i = 0; i < numPoints; ++i) 
    {
        double t = i * dt;
        double originalValue = original(t);
        double reconstructedValue = reconstructSignal(samples, samplingRate, t);
        mse += std::pow(originalValue - reconstructedValue, 2);
    }
    mse /= numPoints;
    return mse;
}

int main() 
{
    double duration = 1.0;  // 信号持续时间
    double highestFrequency = 10.0;  // 最高有效频率

    double validSamplingRate = 2.5 * highestFrequency;
    std::vector<double> validSamples = sampleSignal(originalSignal, validSamplingRate, duration);
    double validMSE = calculateMSE(originalSignal, validSamples, validSamplingRate, duration);

    double invalidSamplingRate = 1.5 * highestFrequency;
    std::vector<double> invalidSamples = sampleSignal(originalSignal, invalidSamplingRate, duration);
    double invalidMSE = calculateMSE(originalSignal, invalidSamples, invalidSamplingRate, duration);

    std::cout << "采样率 " << validSamplingRate << " Hz 时的均方误差: " << validMSE << std::endl;
    std::cout << "采样率 " << invalidSamplingRate << " Hz 时的均方误差: " << invalidMSE << std::endl;

    // 输出采样数据到文件，用于可视化
    std::ofstream validFile("valid_samples.txt");
    std::ofstream invalidFile("invalid_samples.txt");
    double dt = 1.0 / 1000;
    for (double t = 0; t < duration; t += dt) 
    {
        validFile << t << " " << reconstructSignal(validSamples, validSamplingRate, t) << std::endl;
        invalidFile << t << " " << reconstructSignal(invalidSamples, invalidSamplingRate, t) << std::endl;
    }
    validFile.close();
    invalidFile.close();

    return 0;
}    